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The study of nonlinear phenomena has changed dramatically during the last two decades or 
so, as an increasing number of once forbidding problems have become amenable to treatment 
by faster and cheaper computers. From coupled anharmonic oscillators to gravitational 
clustering, from plasma physics to the dynamics of phase transitions, numerical simulations 
are often the only tool to probe the physics of complex nonlinear systems [0. 

Typically, we are interested in investigating the behavior of a particular physical system 
described by ordinary or partial nonlinear differential equations. In the present work, focus 
will be mostly on the latter case, which can be thought of as representing systems with finitely 
or infinitely many coupled degrees of freedom. Apart from very few exceptions, such as kink 
solutions for sine-Gordon or cf)^ models 0, nonlinear partial differential equations have no 
analytical solutions. The situation is even worse if we attempt to model realistic behavior by 
coupling the system to an external environment. This external environment often represents 
a thermal or quantum bath, adding an element of stochasticity to the deterministic evolution 
equations. In order to gain some insight into the role of nonlinearities, perturbation theory 
is frequently used. However, examples ranging from the simple pendulum equation to 
critical phenomena during phase transitions remind us that perturbation theory breaks 
down precisely in the region of parameter space where nonlinear effects become predominant. 

The alternative is to address the problem numerically, solving the equations of interest 
using a computer. In the case of partial differential equations, the problem is set up on 
a lattice which represents a particular choice of discretization procedure. For a function 
of c?- dimensional position and time, /(x, t), satisfying some partial differential equation 
with given initial and boundary conditions, we typically construct a (i-dimensional lattice 
of a given geometry, say cubic or triangular, to represent space at a particular instant, 
and replicate it at (usually regular) intervals to represent time. The continuous function 
may then be discretized following well-prescribed rules by which continuous derivatives are 
approximated by finite ratios of the lattice variables 0. 

The use of a spatial lattice introduces two artificial length scales; the 'macroscopic' size 
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of the lattice in each dimension, L, and the 'microscopic' distance between neighbouring 
lattice points, 5x. These length scales provide bounds on the wavelengths of modes which 
can be represented on the lattice, whilst the total the number of lattice points (for cubic 
lattices being = {L/5xY) is the restricted number of degrees of freedom being integrated 
at each time step. Computational physicists (and computers) spend a considerable amount 
of time trying to get around the limitations that these length scales introduce to numerical 
studies of continuum systems. Occasionally, one or other of these limitations may become 
insignificant due to the particular physical behavior of the system; for example, close to the 
critical point of a second order phase transition the divergence of the characteristic length 
scale of the system means that its bulk properties (and in particular its critical exponents) 
are determined by the long wavelength modes alone, doing away with the need for the high 
spatial resolution given by a small lattice spacing 5x [Q. In general, however, since the 
continuum corresponds to the limit L— >oo,5x— >0,A^— >oo, a better approximation is 
obtained from a larger and finer lattice, leading to the notion of the continuum limit of a 
discrete system. For continuum systems described by continuous functions, such as fluids, 
fields, or deformable bodies, our discrete representation should have a well-defined continuum 
limit, i.e., one that is stable as 5x — (at fixed L). Moreover, we should also demand that 
it is a good continuum limit, in that it matches the original continuum system. As discussed 
below, for systems coupled to external environments, even if the continuum limit can be 
achieved on the lattice it is not always clear how to match the lattice results to a continuum 
theory. These two questions — how to achieve a continuum limit in lattice simulations, and 
how to ensure that it is a good limit, in the sense of matching the appropriate continuum 
theory — are the focus of this work. 

For linear systems, achieving a continuum limit does not usually present any difficulties. 
Typically there is a minimal length-scale in the problem which can be used as a guideline for 
the choice of 5x. For example, when solving the wave equation, it is possible to find a small 
enough 5x and show that the same results are obtained if smaller values are used, provided 
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one makes sure the discretization of time is appropriately chosen so that the evolution is 
stable. 

For nonlinear systems, the situation is more complicated. If we think for a moment in 
terms of a Fourier decomposition of the function /(x, t), the effect of nonlinearities is to 
couple different wavelength modes in a nontrivial way; the dynamics of short wavelength 
modes will influence the dynamics of long wavelength modes and vice-versa. Mechanisms to 
handle this problem are sensitive both to the particular system under study and to which of 
its properties are of interest, often seeming to be more an art than a science. For example, 
if we are solely interested in the dynamics of long wavelength modes with slow relaxation 
time-scales, it may be possible to add extra artificial terms to the evolution equations which 
damp the behavior of faster modes. For situations in which nonlinear fields are coupled 
to an external environment with stochastic properties, say a thermal (or quantum) bath, 
a detailed investigation of how to approach the continuum limit on the lattice is lacking. 
This does not imply that this problem has been completely overlooked, but that it may have 
received less attention than it deserves. 

In the context of classical field theories at finite temperature there has been some work 
on obtaining such a continuum limit. For example, Parisi |^ suggested the addition of 
renormalization counterterms, a proposal then implemented by Alford and Gleiser in 
the context of 2-dimensional nucleation studies (albeit with a somewhat ad hoc match to a 
continuum theory), and by Kajantie et al. [|| in lattice gauge simulations of the electroweak 
phase transition. Alford and Gleiser in particular showed that neglecting lattice spacing 
effects in the numerical determination of nucleation rates can lead to severe errors in the 
measured values. This conclusion is not particular to systems exhibiting metastable states, 
but to any nonlinear field model in contact with external stochastic environments. Thus, 
the issues that are raised here are of concern to a wide range of physical systems modelled 
through the separation of system and environment, from quantum field theories to effective 
field theories describing condensed matter systems. 
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Even if a continuum limit can be achieved on the lattice, we must still ensure that the 
numerical results correspond to the appropriate continuum theory. In general, the coupling 
to a stochastic environment modifies the effective lattice theory, which cannot be naively 
matched to the original continuum model. The question then becomes what theory is the 
lattice simulating, and can we extract it in a self-consistent way? These questions will be 
addressed below in the context of two continuous nonlinear models in 2+1 dimensions, one 
temperature independent and the other temperature dependent (the well-known Ginzburg- 
Landau model). Both models describe phase transitions in the Ising universality class. 
Extensions to d + 1 dimensions should be straightforward. 

Formulating continuum models on a lattice: The issues 

Consider a single scalar field 0(x, t) in a potential Vo{(f)) which may or not be temperature 
dependent. This potential can model interactions of with itself and with other fields. For 
example, a linear term of the form 4>7i is often used to represent the coupling of 4> to an 
external magnetic field for models of ferromagnetic transitions. In this report, focus will 
be on potentials which are simple polynomials of even power in 4>, although our approach 
is equally valid for potentials with odd powers of 4>, typical of nucleation studies. The 
Hamiltonian for this system is, (in units oi c — — i) 

The field 4> can be thought of as representing a scalar order parameter in models of 
phase transitions in the Ising universality class, such as ferromagnets, binary fiuid mixtures, 
metal alloys, or in studies of domain wall formation in cosmology. As such, it is convenient 
to model its dynamics in contact with a heat bath by means of a generalized Langevin 
equation, 

_| = VV-,^-^ + «x,t) , (2) 

where the viscosity coefficient r] is related to the stochastic force of zero mean ^(x, t) by the 
fluctuation-dissipation relation. 



(e(x, t)e(x', t')) = 2r^T<5(x - x')(5(t - t') . (3) 

This approach guarantees that will be driven into equilibrium, although the time-scale 
1]^^ is arbitrary. It has been extensively used in numerical simulations of thermal creation of 
kink-antikink pairs nucleation spinodal decomposition |TI|, and pattern formation 



in the presence of external noise , to mention but a few examples. Note that in the high 
viscosity limit the second-order time derivative can be neglected, as is common practice in 
systems with slower dynamical time-scales. 

The next step is to discretize this system and cast it on a lattice. Using a standard 
second-order staggered leapfrog method we can write, 

= <Pi,m + 8t4>i,m+ll2 (4) 

where z-indices are spatial and m-indices temporal, overdots represent derivatives with re- 
spect to t and primes with respect to 0. The discretised fluctuation-dissipation relation now 
reads 

= 2^^^^ (5) 

so that 



'^'^^-G^,m (6) 



where is taken from a zero-mean unit-variance Gaussian. 

Note that as a first guess we have used Vo(0) in the lattice formulation of the model. 
Is this the correct procedure? It is well-known that classical field theory in more than one 



spatial dimension is ultraviolet divergent, the Rayleigh- Jeans ultraviolet catastrophe p3 
Formulating the theory on a lattice takes care of the problem, as a sharp momentum cutoff is 
introduced by the lattice spacing 5x, with A = n /5x. However, a finite lattice spacing creates 
two difficulties. First, the lattice theory is coarse-grained on the scale 5x\ in other words, 
the lattice theory is not equivalent to the continuum theory we started with, and our results 
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will depend on 6x, unless this dependence is handled by a proper renormalization procedure. 
Second, if the lattice theory is not equivalent to the continuum theory we started with, to 
what continuum theory is it equivalent to? Fortunately, there is a well-defined procedure 
that addresses both difficulties at once. Within its validity, it is possible to establish a 
one-to-one correspondence between lattice simulations and field theories in contact with 
stochastic baths. 

Formulating continuum models on a lattice: The procedure 

In order to recover the continuum limit on the lattice we must eliminate any dependence 
on the cutoff. The coupling to the heat bath will induce ffuctuations on all possible scales. 
Since the cutoff sets the scale for the smallest possible spatial ffuctuations in the system, we 
may incorporate the effects of all ffuctuations down to the smallest scale using perturbation 
theory. Thus, the lattice theory must be equivalent to a continuum theory with a sharp 
ultraviolet cutoff. For classical ffeld theories, the one-loop corrected effective potential with 
a large momentum cutoff is given by , 



= ^0 + /o (27rp^^ counterterms . (7) 

These theories describe ffuctuations with huj <^ ksT. In semi-classical language, the ex- 
citations of the ffeld contain many fundamental quanta. Note that there is a one-to-one 
correspondence between classical statistical ffeld theory in d + 1 dimensions and Euclidean 
quantum ffeld theory in d dimensions. While the loop expansion is in powers of T for the 
former, it is in powers of h for the latter. For d = 2, the only divergences are at one loop, 
although higher loops can generate ffnite terms which modify the effective Hamiltonian. The 
dependence on the cutoff A can be handled by introducing proper counterterms. 
Integration gives. 



ovr 



In 



A2 



-|- counterterms . (8) 



The form of Vq will determine the counterterms needed to cancel the dependence on A. For 
polynomial potentials of order 0", one typically needs counterterms up to order In 



the case of interest here, degenerate double-well potentials, only one quadratic counterterm 
is needed, of form acf?, with a constant. As usual, the value of a is fixed by imposing 
a renormalization condition. Because of the logarithmic divergence, the renormalization 
condition must be imposed at some energy scale M, which is chosen to be. 



M) = 1/o"(0 = VM) . 



(9) 



The renormalized one-loop corrected potential is then. 



ovr 



+ 



T 



+ 



V" 



(10) 



cf>=VM 



The above procedure incorporates thermal fluctuations to the original potential Vo(0) 
at some energy scale M to one-loop order. As with any perturbative approach, it will 
break down wherever large amplitude fluctuations are present, and in particular close to the 
critical point Tc. Although there are techniques to improve the perturbative expansion in 
the neighborhood of the critical point, such as ^-expansion methods (not too reliable 
for 2-d), in this work we will concentrate on the matching of the continuum theory to the 
lattice simulation in regions of the parameter space where the one-loop calculation is valid. 
Close to criticality the theory of Eq. breaks down, and we restrict our investigation to 
the extraction of the critical exponent controlling the divergence of the order parameter. 

How is this continuum theory matched to the lattice simulation? The procedure we 
propose is quite simple. Since the continuum theory above incorporates fluctuations from 
momentum scales up to A, we write the lattice potential as. 



(11) 



where a is fixed by the renormalization condition in the continuum, but with A = 7t/6x. 
That is, 



^att(0) = K) + 



T 



Vr In 
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(12) 



<f>=VM 



As we show below, this procedure takes care of the two problems raised by formulating the 
continuum theory on the lattice, namely, the dependence of lattice results on lattice spacing 
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and the matching of the lattice theory to the continuum at some renormahzation energy 
scale M. The generic emergence of a good continuum limit from Eq. |T2| is the central result 
of this work. 



Applications 

We will apply the above procedure to two cases, with potentials which are tempera- 
ture independent and temperature dependent, respectively. Consider first the temperature- 
independent potential, 

Vo(0) = -^mV' + . (13) 



Choosing the renormahzation point to be (j)RN = y ^^^^^^ , the renormalized continuum 
potential is, from Eq. 

VM = -2-^^ + 4^0^ + ^ + 2 j 0^ - -(3A0^ - m^) In [-^^ J • 

(14) 

It is convenient to introduce dimensionless variables (because there is no h in this theory, m 
has dimensions of (length)"^ while has dimensions of (energy) ^Z^), x = xm, i = tm, = 
(/)A^/^m~^, f] = r]m~^, 9 = T\m~'^, M = Mm~^, A = Am~^. From the discussion in the 
previous section, the lattice-spacing independent lattice potential is, using dimensionless 
variables (and dropping the tildes), 

Viatt(0) = --0' + -0^ + - ( In ( + ] 02 . (15) 



2 4 47r\ Vvr 



Fig. 1 shows the impact of the added counterterm to the lattice results. We display the 
time evolution of the spatially averaged field, (p = \ I dAcp, starting from a broken sym- 
metric phase = — 1, without the counterterm (Fig. la) and with the counterterm (Fig. 
lb). The parameters 6, M, and physical lattice size L, were kept fixed, and only the lattice 
spacing 6x was varied. (Throughout this work we keep the viscosity coefficient = 1 as we 
are only interested in final equilibrium quantities.) Clearly, omitting the counterterm leads 
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to severe lattice spacing dependence of the results, even to the point of having symmetry 
restoration. Experiments varying 9 and M showed that the procedure is robust, with excel- 
lent 5a;- independence being achieved, even close the critical point, as long as the expansion 
parameter 6'/87r ^ 1. 

The next step is to compare the lattice results with the continuum models of Eq. ^ 
in their domain of validity. Being perturbative, we expect the continuum models to break 
down when the fluctuations become large, at high temperatures or close to the critical point. 
By contrast, the lattice models incorporate fluctuations up to the limiting size L, and so 
may remain valid even when the continuum models break down. The continuum potential 
gives a prediction for the critical temperature of 

27r 

' ~ 3(l + M-2 + lnM) ■ 

Note that 9c has its maximum value at = 2; as we move away from this point in 
either direction 9c decreases, and we should expect perturbation theory to continue to be a 
valid approximation closer and closer to the critical point. Ultimately, however, the phase 
transition is nonperturbative, the field fluctuations become large, and perturbation theory 
must fail. Fig. 2 shows the variation in the equilibrium mean field value 0eq with temperature 
6*, squares from the lattice and lines from the continuum, for values of the renormalization 
energy-scale M = 0.1 (Fig. 2a), M = ^2 (Fig. 2b), and M = 10 (Fig. 2c). The 
discontinuities in the continuum are related to the concavity of the corrected potential 
between the inflection points, which gives rise to an imaginary part. As shown by Weinberg 



and Wu [|T6|, the imaginary part of the potential represents unstable physical states typical 
of the process of phase separation; the figure shows only the real part of the corrected 
potential. There is indeed excellent agreement at low temperatures, which is progressively 
lost as the temperature increases. 

At the one-loop level, perturbation theory is equivalent to mean field theory. Close to 
the critical point, where mean field theory breaks down, we expect the equilibrium value of 
to diverge as a power law, 
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0eqOc|— ^1 (17) 

with the critical exponent (3 = 1/2 for mean field theory and /3 = 1/8 for the 2-d Ising model. 
Figure 3 shows the behavior of the lattice and continuum equilibrium mean field values 0eq 
with reduced temperature 9r = {9c — 0)/9c for M = 0.1 — squares being results from 
the lattice simulations, triangles the predicted behavior from the continuum, and the lines 
indicating the two slopes (3 = 1/2 and /3 = 1/8. We see that the continuum perturbation 
theory behaves as a mean field theory, whilst the lattice theory in the neighborhood of the 
critical point is in the universality class of the 2-d Ising model as expected. 

We now consider the case of a temperature-dependent potential. The goal is to show 
that the above procedure works equally well in this case; both lattice-spacing independence 
and the matching to a continuum theory can be achieved in a consistent way. Coupling a 
temperature-dependent potential to a heat bath does not necessarily imply a double counting 
of the thermal degrees of freedom. The choice of potential Vq simply refiects different physical 
models. For example, one may include phenomenological temperature-dependent terms in 
Vq, as in the Ginzburg-Landau model, or may obtain temperature corrections by integrating 
out from the partition function either other fields coupled to or short wavelength modes 
of the field itself |T^. In either case, the heat bath may then be representing stochastic 



forces not included in the integration process, or simply an external environment coupled to 
(j) phenomenologically, which drives the system to its final equilibrium state. As an example, 
we choose the Ginzburg-Landau potential, 

= la{T - Z)<j>' + ^A0^ , (18) 

where the prime is a reminder that the critical temperature has an arbitrary value in the 



mean field model. Fixing the renormalization energy scale at (pRN = \l ^ '^3\ ~^ > 
renormalized continuum potential becomes, 

, , 1 , „ . 1 4 SAT / M^-a(T- T') \ „ 



2 ' ^'"^ 4 87r \ M2 
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-£K + «(r-r^,].n(?^^l±^). (19) 

Following the same steps as before and arbitrarily setting 9'^ = 1, this theory is matched on 
the lattice to 



1 1 "^f) 



(20) 



Fig. 4 compares the lattice results without (Fig. 4a) and with (Fig. 4b) the renormalization 
counterterm. The prescription to obtain lattice-spacing independence works equally well in 
this case. Fig. 5 again compares the lattice simulations (squares) and the continuum model 
(lines) for renormalization scales M = 0.1 (Fig. 5a ), M = V2 (Fig. 5b), and M = 10 
(Fig. 5c). For low temperatures excellent agreement is obtained, as in the temperature 
independent case. Note that this also confirms that our model has not been 'twice-cooked'; 
had it been, no such agreement would be possible. Finally, in Fig. 6, we show the critical 
behavior of the lattice (squares) and continuum (triangles) for M — 0.1. Again the lattice 
obtains the Ising critical exponent, (3 = 1/8, close to criticality. 

In summary, we have presented a self-consistent method to match lattice simulations to 
nonlinear field theories in contact with an external stochastic environment. This approach 
is of potential interest in a wide range of physical problems, from noise-induced pattern- 
forming instabilities and phase separation in condensed matter physics to symmetry breaking 
in high energy physics and cosmology. It was shown that adding the right renormalization 
counterterms to the lattice potential provides a good continuum limit, independent of the 
lattice-spacing and matching the appropriate continuum theory. That this matching breaks 
down at high temperatures and/or close to a critical point is not surprising, as it reflects 
the limitations of perturbation theory in probing critical phenomena quantitatively. The 
procedure was demonstrated to work well for a large class of widely-used potentials — both 
temperature independent and dependent — and over a wide range of the renormalization 
energy scale M. 
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List of Figures 

Figure 1. The time evolution of the mean field 0(i) at five different lattice spacings 
6x — 0.125, 0.25, 0.5, 1.0 and 2.0 for the temperature independent potential — (a) without 
the renormalisation counterterms added {6x increasing downwards), and (b) with the renor- 
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malisation counterterms added. 

Figure 2. The variation in the equihbrium mean field 0eq with the dimensionless temper- 
ature 6 from the lattice (squares) and the continuum (lines) for the temperature independent 
potential — (a) for M = 0.1, (b) for M = v^, and (c) for M = 10. 

Figure 3. The variation in the equilibrium mean field ^eq with the reduced dimensionless 
temperature Or from the lattice (squares) and the continuum (triangles) for the temperature 
independent potential. The dashed lines have slopes of 1/8 and 1/2. 

Figure 4. The time evolution of the mean field 4>{t) at five different lattice spacings 
5x = 0.125, 0.25, 0.5, 1.0 and 2.0 for the temperature dependent potential — (a) without 
the renormalisation counterterms added {5x increasing downwards), and (b) with the renor- 
malisation counterterms added. 

Figure 5. The variation in the equilibrium mean field 0eq with the dimensionless temper- 
ature 6 from the lattice (squares) and the continuum (lines) for the temperature dependent 
potential — (a) for M = 0.1, (b) for M = \/2, and (c) for M = 10. 

Figure 6. The variation in the equilibrium mean field 4>eq with the reduced dimensionless 

temperature 6r from the lattice (squares) and the continuum (triangles) for the temperature 
dependent potential. The dashed lines have slopes of 1/8 and 1/2. 
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